/*
** Two types of polytopes are implemented.
** One is 2D polygon. The other is 3D polyhedron
** Polytope is the underlying geometric structure of tiles/stencil/compenent
*/

#ifndef _POLYTOPE_H_
#define _POLYTOPE_H_

#include <string>
#include <cmath>
#include <sstream>
#include <valarray>

#ifdef WIN32
#include <Windows.h>
#endif
#include <GL/GL.h>

#include <boost/math/constants/constants.hpp>
#include <boost/shared_ptr.hpp>	// for straight skeleton 

#include <CGAL/Exact_predicates_inexact_constructions_kernel.h>
#include <CGAL/Polygon_2.h>
#include <CGAL/Aff_transformation_2.h>
#include <CGAL/Iso_rectangle_2.h>
#include <CGAL/centroid.h>
#include <CGAL/partition_2.h>
#include <CGAL/bounding_box.h>
#include <CGAL/Straight_skeleton_2.h>
#include <CGAL/create_straight_skeleton_2.h>
#include <CGAL/Partition_traits_2.h>
#include <CGAL/Constrained_Delaunay_triangulation_2.h> // polygon sampling

#include "circle_union.h"

namespace Packing_common
{

	struct Parameter
	{
		typedef CGAL::Aff_transformation_2<CGAL::Exact_predicates_inexact_constructions_kernel> Transformation_2;
		double tx;
		double ty;
		double theta;
		double k;
		Parameter(): k(1.0), theta(0.0), tx(0.0), ty(0.0) {}
		Parameter(double _k, double _theta, double _tx, double _ty) : k(_k), theta(_theta), tx(_tx), ty(_ty) {}
		Parameter& operator*=(double f)
		{
			k = 1.0 + (k - 1.0)*f;
			theta *= f;	tx *= f; ty *= f;
			return *this;
		}
		inline Transformation_2 to_transformation() const
		{
			double ksin = k*std::sin(theta), kcos = k*std::cos(theta);
			return Transformation_2(kcos, -ksin, tx, ksin, kcos, ty);
		}
	};
	inline Parameter operator*(const Parameter& p, double f)
	{
		return Parameter(1.0 + (p.k - 1.0)*f, p.theta*f, p.tx*f, p.ty*f);
	}
	/* polytope as base class */
	class Packing_polytope
	{

	public:
		
		// query
		virtual unsigned int dimension() const = 0;
		// unsigned volume
		virtual double volume() const = 0; 
		// perimeter for 2D polygon, surface area for 3D polytope
		virtual double perimeter() const = 0;

		virtual void gl_draw(const GLfloat *fill_color) const = 0;
		virtual void redraw(const GLfloat *fill_color) const = 0;
		virtual void gl_draw_bbox() const = 0; // draw bounding box

		virtual ~Packing_polytope() { }

	};

	typedef CGAL::Exact_predicates_inexact_constructions_kernel K;
	typedef CGAL::Polygon_2<K> CGAL_polygon;

	class Packing_polygon : public Packing_polytope, public CGAL_polygon
	{
		typedef CGAL::Straight_skeleton_2<K> Straight_skeleton_2;
		typedef CGAL::Exact_predicates_inexact_constructions_kernel::Line_2 Line_2;
	public:

		typedef CGAL::Iso_rectangle_2<K> Iso_rectangle_2;
		typedef CGAL_polygon::Point_2 Point;
		typedef K::Vector_2 Vector;
		typedef CGAL::Aff_transformation_2<K> Transformation_2;
		typedef CGAL_polygon::Vertex_const_iterator Vertex_const_iterator;


		friend Packing_polygon transform(const Packing_polygon::Transformation_2& t, const Packing_polygon& p);

	public:
		//Packing_polygon(const std::string& poly_fm_filename);

		template <class ForIter> Packing_polygon::Packing_polygon(ForIter first, ForIter last) : CGAL_polygon(first, last), id(-1)
		{
			//update_skeleton();
			scale_factor = 1.0; 
			large_neighbor_space = false;
		}
		Packing_polygon() : CGAL_polygon(), id(-1) { scale_factor = 1.0; large_neighbor_space = false; }
		~Packing_polygon();

		//query
		unsigned int dimension() const { return 2; }
		void lib_id(int id_no) { id = id_no; }
		int lib_id() const { return id; }
		double volume() const { return std::fabs(CGAL_polygon::area()); }
		double perimeter() const ;
		Iso_rectangle_2 bbox() const;
		inline double get_scale_factor() const { return scale_factor; }
		inline void set_scale_factor(double f) { scale_factor = f; }
		Point centroid() const { return CGAL::centroid(vertices_begin(), vertices_end(), CGAL::Dimension_tag<0>()); }

		inline bool large_neighbor_room() const { return large_neighbor_space; }
		inline void large_neighbor_room(bool val) { large_neighbor_space = val; }

		Vertex_const_iterator vertices_begin() const { return CGAL_polygon::vertices_begin(); }
		Vertex_const_iterator vertices_end() const { return CGAL_polygon::vertices_end(); }
		const Circle_union& medial_circles() const { return skeleton_circles; }

		//draw
		void gl_draw(const GLfloat *fill_color) const;
		void gl_draw_bbox() const;
		void gl_draw_circles(const GLfloat *fill_color) const;
		void redraw(const GLfloat *fill_color) const;
		void update_skeleton(); //{ skeleton = CGAL::create_interior_straight_skeleton_2(vertices_begin(), vertices_end()); }

		// geometry
		double squared_distance(const Point& p) const;
		// return the gradient of squared distance with x and y, also the closest point
		double squared_distance(const Point& p, Vector& gradv, Point& closest_point) const; 

		//translation
		void translate(double dx, double dy);
		void rotate(double angle) ;
		void self_rotate(double angle);
		void scale(double factor); // scale to factor of the original size
		void relative_scale(double factor); // scale to factor of the current size
		void transform(const Parameter& p);
		Packing_polygon& operator*=(const Transformation_2& t);

	private:
		mutable bool update_display_list;
		mutable GLuint display_list_id;

		double scale_factor;

		int id;

		boost::shared_ptr<Straight_skeleton_2> skeleton;
		Circle_union skeleton_circles;

		bool large_neighbor_space;
	};

	class Polygon_sampler
	{
		typedef CGAL::Exact_predicates_inexact_constructions_kernel K;
		typedef CGAL::Triangulation_vertex_base_2<K> Vb;
		typedef CGAL::Constrained_triangulation_face_base_2<K>  Fb;
		typedef CGAL::Triangulation_data_structure_2<Vb, Fb> TDS;
		typedef CGAL::Constrained_Delaunay_triangulation_2<K, TDS, CGAL::Exact_predicates_tag> CDT;
		typedef CGAL::Vector_2<K> Vector_2;
		typedef CGAL::Iso_rectangle_2<K> Bounding_box;

	public:
		typedef CGAL::Point_2<K> Point_2;

	public:
		Polygon_sampler(const Packing_polygon& polygon_) ;
		void inner_sample(size_t ref_pnt_nb, std::vector<Point_2>& pnts);
		void inner_sample(double ref_step, std::vector<Point_2>& pnts) ;
		void edge_sample(size_t ref_pnt_nb, std::vector<Point_2>& pnts);
	
	private:
		const Packing_polygon& polygon;
		CDT cdt;
	};

	// the polygon is expressed as a container of segments
	template <class SegmentIter>
	void edge_sample(SegmentIter first, SegmentIter last, size_t ref_pnt_nb, std::vector<Polygon_sampler::Point_2>& pnts)
	{
		typedef Polygon_sampler::Point_2 Point_2;
		pnts.reserve(ref_pnt_nb);
		std::list<double> edge_len;
		double perimeter = 0.0;
		for (SegmentIter it = first; it != last; ++it)
		{
			double len = std::sqrt(it->squared_length());
			edge_len.push_back(len);
			perimeter += len;
		}
		std::list<double>::const_iterator len_it = edge_len.begin();
		for (SegmentIter it = first; it != last; ++it, ++len_it)
		{
			size_t n = ref_pnt_nb*(*len_it)/perimeter;
			n = std::max<size_t>(n, 1);
			Point_2 src = it->source(), dst = it->target();
			for (size_t i = 0; i < n; i++)
			{
				double w = static_cast<double>(i)/n;
				Point_2 p((1 - w)*src.x() + w*dst.x(), (1 - w)*src.y() + w*dst.y());
				pnts.push_back(p);
			}
		}
	}

	Packing_polygon transform(const Packing_polygon::Transformation_2& t, const Packing_polygon& p);
	// input function template for Packing_polygon set from files
	template <class InsertIter>
	InsertIter read_polygon_set(const std::string& poly_fm_filename, InsertIter res, const std::string fmt = "poly")
	{
		std::ifstream fin(poly_fm_filename.c_str());
		std::string error_info("Bad input file format :" + poly_fm_filename);
		if ( !fin.is_open() )
			throw Bad_input("Polygon file cannot be openned.\n");
		std::cout<<"Open file: "<<poly_fm_filename<<" as polygon input file..."<<std::endl;
		// backup for memory safety
		typename InsertIter::container_type okay_elements;	
		InsertIter backup_iter(okay_elements);
		try
		{	
			if (fmt == "poly")
			{
				// gross number of polygons
				unsigned int nb_polygons;
				if (!(fin>>nb_polygons))	throw Bad_input(error_info);

				for (unsigned int i = 0; i < nb_polygons; i++)
				{
					unsigned int nb_vtx;
					if ( !(fin >> nb_vtx) )	throw Bad_input(error_info);
					Packing_polygon *np = new Packing_polygon;
					*backup_iter = *res = np;
					for (unsigned int j = 0; j < nb_vtx; j++)
					{
						double x, y;
						if ( !(fin >> x) )	throw Bad_input(error_info);
						if ( !(fin >> y) )	throw Bad_input(error_info);
						np->push_back(Packing_polygon::Point(x,y));
					}
					//np->update_skeleton();
					++res;
					++backup_iter;
				}
			}
			else if (fmt == "txt") // lack of care on exception detection and handling
			{
				int polygon_id = 0;
				while (fin)
				{
					std::string token;
					std::istringstream iss;

					if (!(fin >> token)) break;
					while (token != "PIECE")	
						if ( !(fin >> token) ) // end of file
							break; 

					std::getline(fin, token); // "piece 1"
					std::getline(fin, token); // "quantity"
					unsigned int nb_polygons, nb_vtx;

					//if ( !(fin>>nb_polygons) )	throw Bad_input(error_info);
					std::getline(fin, token);
					iss.str(token);
					iss >> nb_polygons;

					std::getline(fin, token);	// "number of vertices"

					//if ( !(fin>>nb_vtx) )		throw Bad_input(error_info);
					std::getline(fin, token);
					iss.str(token);
					iss.seekg(std::ios_base::beg);
					iss >> nb_vtx;

					std::getline(fin, token); // "vertices (x, y)"
					Packing_polygon **np = new Packing_polygon*[nb_polygons];
					for (unsigned int i = 0; i < nb_polygons; i++)
					{
						np[i] = new Packing_polygon;
						np[i]->lib_id(polygon_id);
					}
					
					for (unsigned int i = 0; i < nb_vtx; i++)
					{
						double x, y;
						if ( !(fin >> x) )	throw Bad_input(error_info);
						if ( !(fin >> y) )	throw Bad_input(error_info);
						for (unsigned int j = 0; j < nb_polygons; j++)
							np[j]->push_back(Packing_polygon::Point(x,y));
					}
					for (unsigned int i = 0; i < nb_polygons; i++, ++res, ++backup_iter)
					{
						*backup_iter = *res = np[i];
						//np[i]->update_skeleton();
					}
					polygon_id++;
				}
			}
			
		}
		catch (Bad_input&)
		{
			// release allocated elements
			typename InsertIter::container_type::iterator it;
			for (it = okay_elements.begin(); it != okay_elements.end(); ++it)
				delete *it;
			throw; //rethrow for external handle
		}
		std::cout<<"Finished reading polygons...\n";
		return res;
	}
}

#endif